function [update] = estimate_symmetric_update_field(currentScale)
% ESTIMATE_UPDATE_FIELD Estimates the symmetric update field for current scale and iteration
%
% [update] = estimate_update_field(currentScale)
%
% Data used is contained in the data struct registration.
% See also REGISTRATION_EXECUTE
%
% INPUT ARGUMENTS
% currentScale                  - Current scale
%
% OPTIONAL INPUT ARGUMENTS
% N/A
%
% OUTPUT ARGUMENTS
% update
%  displacementUpdateForward    - Estimated update field from fixed to
%                                 moving deformed
%  displacementUpdateBackward   - Estimated update field from moving to
%                                 fixed deformed
% (only if transformation model is set to translation or affine)
%  transformationMatrixForward  - Estimate transformation matrix from fixed 
%                                 to moving deformed
%  transformationMatrixBackward - Estimate transformation matrix from moving 
%                                 to fixed deformed

% Copyright (c) 2012  Daniel Forsberg
% danne.forsberg@outlook.com
%
% This program is free software: you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation, either version 3 of the License, or
% (at your option) any later version.
%
% This program is distributed in the hope that it will be useful,
% but WITHOUT ANY WARRANTY; without even the implied warranty of
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
% GNU General Public License for more details.
%
% You should have received a copy of the GNU General Public License
% along with this program.  If not, see <http://www.gnu.org/licenses/>.

global registration

if strcmp(registration.transformationModel,'non-rigid')
    movingDeformed = deformation(registration.moving{currentScale + 1},...
        field_exponentiation(registration.logDisplacementField),...
        registration.interpolation);
    fixedDeformed = deformation(registration.fixed{currentScale + 1},...
        field_exponentiation(backward(registration.logDisplacementField)),...
        registration.interpolation);
    movingDeformedCertainty = deformation(...
        registration.movingCertainty{currentScale + 1},...
        field_exponentiation(registration.logDisplacementField),...
        registration.interpolation);
    fixedDeformedCertainty = deformation(...
        registration.fixedCertainty{currentScale + 1},...
        field_exponentiation(backward(registration.logDisplacementField)),...
        registration.interpolation);
else
    movingDeformed = affine_transformation(registration.moving{currentScale + 1},...
        registration.transformationMatrix);
    fixedDeformed = affine_transformation(registration.fixed{currentScale + 1},...
        inv(registration.transformationMatrix));
    if ~isempty(registration.movingCertainty{currentScale + 1})
        movingDeformedCertainty = affine_transformation(...
            registration.movingCertainty{currentScale + 1},...
            registration.transformationMatrix);
    else
        movingDeformedCertainty = [];
    end
    if ~isempty(registration.fixedCertainty{currentScale + 1})
        fixedDeformedCertainty = affine_transformation(...
            registration.fixedCertainty{currentScale + 1},...
            inv(registration.transformationMatrix));
    else
        fixedDeformedCertainty = [];
    end
end

switch registration.method
    case 'phase'
        [update] = symmetric_phase_registration(...
            registration.moving{currentScale + 1}, ...
            registration.fixed{currentScale + 1}, ...
            movingDeformed, fixedDeformed, ...
            registration.movingCertainty{currentScale + 1}, ...
            registration.fixedCertainty{currentScale + 1}, ...
            movingDeformedCertainty, fixedDeformedCertainty, ...
            registration.transformationModel, ...
            currentScale);
    case 'optical-flow'
        if registration.multiModal
            [update] = symmetric_optical_flow_registration(...
                registration.moving{currentScale + 1}, ...
                registration.fixed{currentScale + 1}, ...
                movingDeformed, fixedDeformed, ...
                registration.transformationModel,...
                registration.multiModal,...
                'numberOfChannels',registration.numberOfChannels(currentScale + 1));
        else
            [update] = symmetric_optical_flow_registration(...
                registration.moving{currentScale + 1}, ...
                registration.fixed{currentScale + 1}, ...
                movingDeformed, fixedDeformed, ...
                registration.transformationModel,...
                registration.multiModal);
        end
    case 'polynomial-expansion'
        if registration.multiModal
            [update] = ...
                symmetric_polynomial_expansion_registration(...
                registration.moving{currentScale + 1}, ...
                registration.fixed{currentScale + 1}, ...
                movingDeformed, fixedDeformed, ...
                registration.multiModal,...
                registration.transformationModel,...
                'numberOfChannels',registration.numberOfChannels(currentScale + 1),...
                'signalModel',registration.signalModel); 
        else
            [update] = ...
                symmetric_polynomial_expansion_registration(...
                registration.moving{currentScale + 1}, ...
                registration.fixed{currentScale + 1}, ...
                movingDeformed, fixedDeformed, ...
                registration.multiModal,...
                registration.transformationModel,...
                'signalModel',registration.signalModel);
        end
end